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Abstract 

The nonnegativity of the determinant of the partial transpose of a two-qubit (4 x 4) density matrix 
(p) is both a necessary and sufficient condition for the separability of p. While the determinant 

-t— > ■ 

of p itself is restricted to the interval [0, ggg], the determinant of the partial transpose (|/5 PT |) can 



range over \—js, off], with negative values corresponding to entangled states. We report here the 



exact values of the first nine moments of the probability distribution of \p | over this interval, 
with respect to the Hilbert- Schmidt (metric volume element) measure on the nine-dimensional 

OS : 

convex set of real two-qubit density matrices. Rational functions C2j(m), yielding the coefficients 



of the 2j-th power of even polynomials occurring at intermediate steps in our derivation of the 



CO 

co i 

m-th moment, emerge. These functions possess poles at finite series of consecutive half-integers 
(m = — |, — ^, . . . , 2j 2 ~ 1 ), and certain (trivial) roots at finite series of consecutive natural numbers 
(m = 0, 1, . . .). Additionally, the (nontrivial) dominant roots of C2j(m) approach the same half- 



X 

^ . niLcgci vaiuco yiu — — — ,— — , 

found-when employed in the one-sided Chebyshev inequality-give an upper bound of 3 " ' ' 



integer values (m = 2 \ - , 2 \ 3 , . . .), as j increases. The first two moments (mean and variance) 



34749 

0.874759 on the separability probability of real two-qubit density matrices. We are able to report 
general formulas for the m-th moment of the Hilbert-Schmidt probability distribution of \p\ over 
[0, 2gg], in the real, complex and quaternionic two-qubit cases. 
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I. INTRODUCTION 



One interesting, and seemingly not immediately obvious consequence of certain well- 



known results of Peres and the Horodecki clan 



is that one only needs to evaluate the 



sign of the determinant of the partial transpose of a two-qubit (4 x 4) density matrix (p) 
to assess the separability of p pHa], rather than checking individually the signs of its four 
eigenvalues (since no more than one eigenvalue of \p PT \ can be negative). If one assigns a 
measure-we will here use the volume element of the Hilbert-Schmidt (Euclidean/flat) metric 

n n 

[7| (cf. [8[)-to the two-qubit density matrices, then, from the associated probability distri- 
bution over the determinant of the partial transpose, one should-in principle-be able to 



derive t 



(cf. 



.0 



re specific and long-sought probability that a two-qubit density matrix is separable 



llj). To attempt to fully characterize such a probability distribution of interest, 



we begin by computing its first several moments (sec. (TT])- (It has been conjectured that 



"most of the information defining a compact 
is usually contained in its first few moments" 



y supported [probability distribution function] 
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14j.) As a complementary exercise, we simi- 



larly analyze-but with considerably less severe computational demands-the Hilbert-Schmidt 
probability distribution over the determinant \p\ itself (sec. Ill ip . The results obtained allow 
us to construct a general formula (1551) for the m-th moment of this distribution. 

We will proceed in the framework of the Bloore (or correlation coefficient) parameteri- 



zation of the 4x4 density matrices is|-17( which allows us (in the generic real two-qubit 
case of immediate interest here) to work primarily in seven dimensions, rather than the nine 
naively expected. Also, in our computations, we will further reparametrize three of the six 
correlations 



1 < z < j < 4, 



Zij G 



-1,1] 



(1) 



in terms of partial correlations [la], allowing certain requisite integrations to be performed 
simply over six-dimensional hypercubes, rather than more complicated domains. (We had 



181 ] to define the 



alternatively attempted to utilize the cylindrical algebraic decomposition 
integration limits (as indicated in ItI, sec. II]) that specify the domain of feasible density 
matrices, directly within the Bloore-type framework, without transforming to such partial 
correlations. However, we encountered certain apparently inconsistent /puzzling results ob- 
tained using Mathematica in this regard.) 
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II. HILBERT-SCHMIDT PROBABILITY DISTRIBUTION OVER \p PT 



The computations of the m-th moment proceeds in two stages. In the first, we perform an 
integration over the six- dimensional hypercube [— 1, l] 6 of the m-th power of a (transformed) 
polynomial (P)-proportional to |p PT |-in seven variables ( 17j, eq. (7)]). (The proportionality 
factor is (p22P33) 2m -) The free (unintegrated) variable is of the form 



PnP44 
P22P33 



where the Pa s are the diagonal entries of p. (In the related study 



as a variable, and in 



0, 



(2) 



v = /i 2 was used 



19j. £ = log p.) We have that (before the transformation to partial 



correlations, yielding P) 



P = -z^fj, 4 + 2z a (zijZ ik + zpzu) p 3 + 22 jk (z^i + z^z^) /i - zf k + (3) 



^ 2 {{ Z kl ~~ l) Z ij ~ 2 (^il^jk + Z ik Zj\) Z k \Zij + z\a 



^jk Z ]\ Z kl -^ik^il^jk^jl -r ^ ik ^ji 

The transformation of the three correlations z$, z^, Zji to partial correlations (z^^, Zji^, zujk) 
takes the form 16] 



- 2z ik z ii z ik z il + zf k (z£ - l) + l) . 



z n -> ZijZ^zu + J zA - IJzL - lz iKj z u + z^JzL - IJz^ - lzji, fc + 



(4) 



Zl - Ulz{\ - Ulzl , - 1 Jz? - 1 - - -- -L . /*2 - 1 ~- . 



1^14,23 + \I Z { j — l^jkA/ % — l^ik,j£jl,A!, 



Zik -> 2yz jk + y^ij - ly 4 - l^ikj, z fl z jk z k i + ^4 - 1^/4 _ l2 .. i>fe 
The jacobian for this transformation is 



J(zy, Z ifc , ZfcZ, Zjfcj, Z iJ>fc ) = {zl - l) (4 - l) (4 - l) yjz^j - lyfi 



2 



(5) 



For the m-th moment (Moment m = the indicated six-dimensional integration of P in 
now reparameterized form P over the hypercube defined by z^j G [—1, 1], Zj k G [—1, 1], Zu G 
[— 1, 1], Zik,j G [—1, 1], Zji k G [—1, 1], ^iijfe £ [ — 1) 1] yields-including a normalization factor of 
3^a— the ("intermediate function") result 



4i(a0 



27 
32^2 



1 /•! /-I /-I /"l /-I 



1 J_i J_i J_i J_l J_i 



(6) 



i7 (^jj, Zjfc, 2;^;, Zjfcj', Zji t k) [P(z{j 1 Zjk-, Zkl? Zikjt Zji^ki Ziijk)] d.Zijdzjkdzkidz-ikjd.Zji^dzujk- 



3 



For the first (m = 1) moment, we have the result (Fig. [T]) 



/l0 ,) = + (7) 
uw 5 125 5' y J 



for the second moment (m = 2) 



20898// 12/x 2 3 
2l/ij 35 875 42875 875 35' U 

and for the third (m = 3), 

_^ 12 _ 5 V° _ 27873// 8 _ 466876/i 6 _ 27873/z 4 _ 54/i 2 _ 1 
3 ^ " ~ln ~ "875 42875 1157625 ~ 42875 ~ "875" _ 21' ^ 

At this point, we omit terms of lower order 2j than 2m, the coefficients of which match the 
coefficients C im _2j{. m )- Then, 

/x 16 584/x 14 278884/x 12 8984// 10 65788454// 
4 ^ ~ 33" + 5775 + 282975 + 4851 + 20543985 + " ' ' ^ ' 

3/t 20 _18/u 18 70881// 6 _2178728// 4 _59472398/i 12 _ 4103383444// 
5 ^ ~~ ~ 143 143 49049 441441 4855851 273546273 + " ' ' ^ ^ 

/„(,,) = 2556^ 545V» 335937V* (12) 

VP; 65 17875 2695 315315 v ' 

, 3273117// 6 , 597414184// 4 + 173821048732// 2 + 



86515 7872865 1771394625 

f _// 8 _ 4298// 6 _ 826637/i 24 _ 165865636// 2 71226035//°_ 

7 ^ " ~"85" ~ 27625 303875 8204625 722007 ^ ' 

1947049760374// 8 93373201818911/i 16 33225665966177656// 4 
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6711055065 167776376625 48487372844625 

T/s 3u 32 6672u 30 12986136// 8 4250871568// 6 3319251741068// 4 

Un) = J— -| ^_ -| ^_ -| ^_ -| C 1- (14) 

323 40375 3674125 121246125 14670781125 v 1 

755365923834768/Z 2 2024301386770232//° 61510285844520752// 8 3853435310162220966/i 
+ 826454003375 + 826454003375 + 14049718057375 + 724564031244625 

and 

_ // 6 _ 9774// 4 _ 651051/i 32 _ 8355664//° _ 18384996780// 8 _ 4848288282648/Z 6 _ 
9 ^ " ~133 56525 145775 146965 39122083 1944597655 

(15) 

133915228926036// 4 61222919937476688// 2 396008663496240078//° 2103161056387491292// 8 

15026436425 2809943611475 10677785723605 47564681859695 

For the nine cases (m = 1,...,9) we have been able to explicitly compute so far, the 
coefficients of the corresponding 4m-degree even polynomials / m (/i), as already indicated, 
are symmetric-for reasons not immediately apparent to us-around the // m term. 
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FIG. 1: The six functions I m (/u),m = 1,...,6. The curves for even m curve upward, for odd m 
downward, with the steepness of the curves increasing with m. 



The constant terms (as well as the coefficients of the fj, term) are expressible as 

3(_]Y m 

Co (to) = C Am (m) = — 

4 (m +2) (to + |) 

Additionally, the coefficients of the second and (4m — 2) terms are 

3(-l) m m(2m(4m - 5) - 15) 



(16) 



C 2 (m) = C Am - 2 {m) 



100 (to-|) (to + |) (to+|)' 
Further, the coefficients of the fourth and (4m — 4) terms are 

3(-l) m TO(2m(2m(2m(8m(6m - 7) + 155) - 13) - 1017) - 315) 



(17) 



C 4 (to) = C 4m _ 4 (m) 



19600 (m - f) (to - \) (to + \) (to + §) 



(18) 

(These results were obtained using the "rate", guessing program of C. Krattenthaler, based 
on Mathematica programming of M. Trott.) Still further, M. Trott was able to obtain the 
result (using the FindSequenceFunction command of Mathematica) 

C 6 (m) = C 4m - 6 (m) = (19) 

(-1)™(to - 1)to(4to(2to(2to(to(4to(20to(4to - 11) + 173) - 4303) + 4733) + 14911) - 9165) - 4725) 

529200 (to - I) (to - f) (m - |) (to + |) (m + §) ' 
From the formulas for these coefficients, it is clear that the numerator of the coefficient 
(C 2 j(m) = C 4m _2j(TO)) of /i 2 - 7 is a polynomial of degree 3j, and the denominator is a poly- 
nomial of degree j + 2. (The denominators are very simple in structure ( l24"l) -as evidenced 
above.). For j = 0, the roots are — | and — |, and as j increases by 1, an additional root 
1 larger in value than the previous smallest is added. Thus, poles occur at the coefficient 
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functions at such half-integers.) Utilizing this observation, we were then able-using simple 
fitting procedures-to move on to obtaining the coefficients C%{m) = C^ m s(m), Cio(m) = 
C4m-io(rn) C 12 {m) = C 4m _i 2 (m), C 14 (m) = C 4m _i 4 (m) and C w (m) = C 4m _i 6 (m)-but not 
yet higher. In studying the roots of these functions, we have detected one quite interesting 
feature. That is, as j increases, the dominant roots of (^ (m) show very strong evidence of 
converging to j — |, the subdominant roots to j — |,...For instance, for j = 8, the dominant 
roots of C 16 m = C 4m _i 6 (m) are 7.49999796, 6.4999352, 5.4980028, 4.4493216, while for j = 7, 
they are 6.5000204,5.500556,4.515944. Such roots would then come increasingly close to 
canceling the near-to- matching poles in the denominators in C 2 j (m) as j increases. 

In the second stage of our procedure to compute the m-th moment, we reverse the substi- 
tution §2§ in these 4m-degree polynomials, multiply the result by the necessarily nonnegative 
factor (p22P33) 2m (the factor (P22P33) had been removed in forming the polynomial P in seven 
variables, proportional to |p PT |) and also by the jacobian corresponding to the transforma- 
tion to Bloore (correlation) variables for the real two-qubit density matrices 
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jac = (P11P22P33P44) 2 • (20) 

The result is, then, integrated over the unit three-dimensional simplex, 

P11 + P22 + P33 + P44 = 1, Pa > 0, i = 1, . . . , 4 (21) 

to obtain the m-th moment. In other words (taking into account the appropriate normal- 
ization factor), and setting p 44 = 1 — p n + p 2 2 + P33, 

1146880 f 1 f x - pu [ x -p^-p™ 

Moment m = ( m = — / / / (22) 

77 Jo Jo Jo 



(P22 P33 ) 2m ( Pi 1 P22 P33 P44 ) 2 I m ( \ PllPU ) dp 33 dp 22 dpi i . 

V A22P33 

We are, in fact, able to perform the indicated symbolic integration, obtaining thereby 

M — ™ - <• - ,™ w f ^- T (tt) 2 r ("5 + 2m + §) 2 Ci[m) (23) 

2293760 _ 9m „ 9 ^fi + 5\ 2 ^f i „ 5 X 2 



7r 2 r(4m+10) W4 - V 2 J V 2 2 



1146880 „ / 5 N 4 
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where the Cj(m)'s are our previously-indicated rational functions ( (jT6|) -( TT9|) ). symmetric 
about 2m. These (rational functions) Cj(m)'s themselves are the ratios of polynomials in m 
of degree y divided by the term (using the Pochhammer symbol, as well as rising factorials 
for gamma functions with half-integer arguments) 

(l-i \ r(m + |) 25+ 2 (2m + 3)!! , lS 

denominator^ (m)) = V m] = — t — — rr = '— (24) 

K K JJ V 2 J 1+2 r(-f + m+i) (-i + 2m-l)!! 1 ; 



1 - k, 



= n' fc= - 2 ,o,.( m + . ; 

For i = 4, by way of example, this gives us ffTB"]) 

m ~£) ( m_ i) ( m+ ( m+ 0- (25) 

On the other hand, the numerators of the Cj(m)'s for m > have zero as a trivial root, and 
for m > An, trivial roots 0, . . . n. 

Again, converting gamma functions with half-integer arguments to rising factorials, we 
have, equivalent ly, that 

Momentm = C m = (26) 

35 r[4m + 1Q] [((2m + 3)!!) C 2m (m) + 2S^ 2 4 ... ((3 + *)<<(3 - i + 4m)!!) ( 

Pursuant to these formulas, the first moment (mean) of the Hilbert-Schmidt probability 
distribution of |p PT | over the interval [— ^g] is (departing from the convention of denoting 
moments by /i, since that symbol has been employed in our earlier studies 13, ll9] and above 
©) 

Ci = — = w -0.0011655, (27) 

S1 858 2-3- 11 -13 v ; 

falling within the [negative] region of entanglement. Then, successively, the ([necessarily] 
decreasing in absolute value) raw (non-central) moments are 

27 3 3 

Co = = -tt; w 0.0000108462, (28) 

S2 2489344 2 10 -11-13-17 ' V ; 

8363 8363 7 , , 

G = = — - — w -1.2629773 ■ 10~ 7 , 29 

53 66216550400 2 13 • 5 2 • 7 ■ 11 ■ 13 • 17 ■ 19 ' V ; 

C' = 21859 = 21859 „ 2 09311 • 10- 9 (30) 

54 10443295948800 2 17 • 3 - 5 2 - 11 - 13 - 17 - 19 • 23 ' ' v ; 

23071 23071 n 

-4.27531 • 10 -11 , (31) 



539633583390720 2 18 ■ 3 ■ 5 • 7 2 • 13 ■ 17 ■ 19 ■ 23 ■ 29 
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3317321 7 -43 -103 -107 , 



3253917653076541440 2 28 • 3 • 5 • ll 2 • 17 • 19 • 23 • 29 • 31 



1.01949- 10' 12 , (32) 



c > = 419856257 _ 

7 15366774022001834065920 1 ' 

43 ' 2179 ■ 4481 --2.73223- 10-", 



2 30 • 3 4 • 5 • 11 • 13 • 17 • 19 • 23 • 29 • 31 • 37 



c; = = (34) 

S8 21117403549591928832000 v ; 

109 ■ 155461 _ aQ2431 _ 



2 33 - 3 - 5 3 -ll-19-23-29-31-37-41 
and (requiring four days of computation on a MacMini machine) 



6102620963 

^ 9 ~ ~ 240565904621616585139814400 ~ ^ ' 

19 • 199 • 1614023 17 

~ -2.53678 • 10~ 17 . 



2 37 • 3 • 5 2 • ll 3 • 13 • 23 • 29 • 31 • 37 • 41 • 43 
(After four weeks of uninterrupted computation, we did not succeed, however, in determining 

CioO 

Interestingly, the sequence of denominators immediately above (in apparent contrast to 
that of the numerators) appears to be "nice" in that the number of their prime factors do 
not grow rapidly, but rather linearly. This is a strong indication of the possible existence of 



a "closed 
factorials 



orm" , that is an expression which is built by forming products and quotients of 
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fn. 12]. 

The skewness (71) of the Hilbert-Schmidt probability distribution over |p PT | is negative 
(as well as all moments listed of odd order), that is, -3.13228-so, the left tail of the distribu- 
tion is more pronounced than the right tail-while its kurtosis (72), a measure of " peakedness" 
is quite high, 17.6316. (Higher kurtosis indicates that more of the variance is the result of 
infrequent extreme deviations than frequent modestly sized deviations.) From the first two 
moments, we obtain the variance 

2 = 30397 ^ 9 487838 . 10 -e (36) 
3203785728 v ; 

Application of the standard- form one-sided Chebyshev inequality {22] (we perform a 
linear transformation, so that negative values of |p PT | are mapped to [0,1]), then, yields an 
upper bound on the Hilbert-Schmidt separability probability of the real two-qubit density 
matrices of fjfff « 0.874759. This is a substantially weaker upper bound than that of 
IjII ~ 0.537619 established in [19], by enforcing the nonnegativity of pairs of 3 x 3 principal 



minors of the partial transpose, as well as weaker than ps 0.76854, obtained by requiring 



19] . 



the nonnegativity of all six 2x2 principal minors of the partial transpose 

Using the general formulae for the coefficients ( (116j) - (118p ) -derived above, using the "rate" 
program of C. Krattenthaler-of the zero-th, second, fourth (and by symmetry) the 4m, 4m— 2 
and 4m — 4 powers of \x in the intermediate functions I m (f/,), we have been able to obtain 
the exact contribution of the associated six terms to the m-th moment. This contribution 
is the product of the two factors 



945(-l) m (m(2m(2m(2m(40m(6m - 5) - 169) + 101) - 495) - 9) + 27) 

and 



2 (16m 4 -40m 2 + 9) 
256r(2m+|) 2 2~ 8 "T(4m + 



(37) 



+ 



7rr(4m + 10) (m + 2) (4m + l) 2 (4m + 3) 2 (4m + 5) 2 (4m + 7) 2 (4m + 9)r(2m + 4) 2 

(38) 

For m = 3, the ratio of the true/known moment to the product of these two factors is 1.05766, 
increasing monotonically, in a quasi-linear manner, to 1.94638 for m = 9. Extending this 
form of analysis/approximation to take into account the exact formulas we have also so far 
obtained for Cijim) = C , 4 m _2j(m), for j = 3, ... ,8, we can reduce this ratio from 1.94638 

t 170368623463798669312 ^ i QQE.14 
164584930558733068259 ~ 



III. HILBERT-SCHMIDT PROBABILITY DISTRIBUTION OVER \p\ 

The determinant of a 4 x 4 density matrix itself is restricted to a smaller range [0, ^] 
than that-[— ^, ^gj-of its partial transpose. We have computed the initial moments of the 
Hilbert-Schmidt probability distribution of \p\ over this interval, where p corresponds to a 
generic real two-qubit system. (Doing so involves only a series of three-fold integrations [?], 



23j -three being the number of independent eigenvalues of a 4 x 4 density matrix-rather than 
the nme-fold [6+3] integrations used above for the moments of the probability distribution 
of |p PT |.) The first moments (q) are 

g = = (2 4 ■ 11 ■ 13)" 1 w 0.000437063, (39) 
2288 

' (2 10 • 11 ■ 13 • 17)" 1 w 4.01712- 10~ 7 , (40) 



2489344 
1 



1891901440 

9 



(2 13 • 5 • 11 • 13 • 17 • 19)" 1 ps 5.28569 • lO" 10 , (41) 



and 



C'a = = (2 17 • 5 2 • 11 • 13 • 17 • 19 • 23) _1 « 8.61797 • 10~ 13 , (42) 

S4 3481098649600 V ; V J 

(I = = (2 21 • 3 • 5 • 7 • 13 • 17 • 19 • 23 • 29)" 1 « 1.62147 • 10" 15 , (43) 

S5 616724095303680 V ; V 1 

C« = « 3.38054 • 10~ 18 , (44) 

S6 295810695734231040 ' v ' 

(I = « 7.61383 • KT 21 , (45) 

S7 131339948905998581760 ' v ; 

(' = « 1.82132 • 10~ 23 , (46) 

S8 54905249228939014963200 v ; 

Co = « 4.57255 • 10" 26 . (47) 

S9 21869627692874235012710400 v 7 

Cin = « 1.19439 ■ 10~ 28 (48) 

S1 ° 8372488740021088953229639680 v ; 

?' — ~ 3 22484 • 10 -31 (49) 

su 3100931560849549878551107338240 ' v ' 

(' = « 8.95482 • 10" 34 (50) 

S12 1116717015021019439340374162669568 v ; 

3 

— ~ 2 54723 • 10~ 36 (51) 

Sl3 1177747849688102259981247358666014720 ' v ; 

?' — ~ 7 39882 • 10~ 39 (52) 

S14 135156718750413942110951421022086103040 ' v ; 

P ^ 2 18883 10~ 4 ^ f53) 

~ 45686458852962503761039927438910834081792 ~ ' ' 1 ' 



IV. RELATIONS BETWEEN THE TWO SETS OF MOMENTS 

Similarly to the sequence (( I27p - (j35|) ) of denominators of the earlier set of moments pre- 
sented, the prime factors of the denominators in this latter set of moments (and nine more 
we have been able to compute) do not grow rapidly, indicative of the possibility of a closed 
form for them. (Of course, the numerators here-mostly l's with two 3's-unlike the earlier 
ones, are comparatively well-behaved.) For instance, 

C24 = (2 98 -3 2 -5 2 -ll-13 2 -17-19-29-31-53-59-61-67-71-73-79-83-89-97-101-103)" 1 « 6.66035-10" 64 . 

(54) 

Interestingly-but for the cases m—1, 5-the powers to which 2 is raised in the denomina- 
tors of the m-th entries (m = 1, 9) of the two moment sequences ( ff271) - fl35|) and fl39l) - fl53l) ) 
coincide. (In our two sequences, none of the powers of 2 occurring in the denominators is a 
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FIG. 2: Fit of a twenty-four degree polynomial to the first twenty-four moments of the Hilbert- 
Schmidt probability distribution over |p|, where p is a generic real two-qubit density matrix 

square number.) More strikingly still, taking 11 to be the fifth prime number, the highest 
prime occurring in the denominator of the m-th member of each of the two sequences is the 
same (5+m)-th prime. 

In Fig. |2l we display the fit of a power series in \p\ of degree twenty- four (with twenty- five 
unknown coefficients) to the computed first twenty-four moments. (We also-giving us the 
needed twenty-fifth constraint-require the "zeroth" moment to be 1, as mandated for any 
probability distribution.) No nonnegativity constraints were, however, imposed and certain 
slight incursions into negative regions result. (The "probability" mass below the \p\ axis 
is 0.000397, while above it, the mass is 1.000397.) The distribution is clearly peaked at 
\p\ = 0, the locus of the degenerate (pure, pseudo-pure,...) states, those having at least 
one eigenvalue zero. In Fig. [3j we attempt an alternative reconstruction of this probability 



24 eq. (6)], giving us a 



distribution using the stable approximation method advanced in 
sequence of plateaus. 

Further, inputting the first twenty-four moments computed into the FindSequenceFunc- 
tion command of Mathematica, we obtained the formula (again employing the Pochhammer 
symbol) 

2 1 - 8 -(l) m (j) m 945v^2- 8m - 4 r(2m + 2) 
Cm (-+2)(^) w (fL (m + 2)T(2 m + ±) ■ 
The associated moment- generating function is the generalized hypergeometric function 

M(t) = 4032 (3 F 3(|,l,l;|,2,|;^g) - 1) ^ 

The inverse Fourier transform of the associated characteristic function, that is M(it), then, 

11 
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FIG. 3: Reconstruction of the Hilbert-Schmidt probability distribution over \p\, based on its first 



twenty-four moments, using the stable approximation method advanced in 



24j, eq. (6)] 



should yield the Hilbert-Schmidt probability distribution over \p\ for generic real two-qubit 
density matrices p. (However, we have not been able to explicitly evaluate it.) 

In Fig. HI we display the fit of a power series in \p\ PT of degree nine to the computed first 
nine moments (( 127]) - (I35]) ) (cf. Figs. 1, 2]) of the Hilbert-Schmidt probability distribution 
over |p PT |, where p is a generic real two-qubit density matrix. No nonnegativity constraints 
were, however, imposed and considerable incursions into negative regions now result. (Such 
negativity can be obviated through the use of maximum-entropy, spline-fitting and other 
methodologies 14], |24j, and we intend to explore such directions. The use of the stable 



approximation approach 



24| used to produce Fig. [3] was not insightful with our relatively 



small number of explicit moments.) 

Since the plotted distribution (Fig. H]) appears to be unimodal, one can presumably use the 
computations of the first and second moments above to isolate the mode of the distribution 
within the interval eq. (13)] 



30397 
51 



30397 



51 



858 



4576 ' 4576 



858 



{-0.00650062, 0.00416962}, 



(57) 



containing |p PT | = 0. Narrower intervals containing the mode can be obtained using higher- 



order moments and the associated Hankel determinants 



m, Thm. 3.2]. 
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-0.0+0 — =03508 -0.006 -0.004 -0.002 



0.002 \ 0.004 



IP I 



FIG. 4: Fit of a nine-degree polynomial to the first nine moments of the Hilbert-Schmidt proba- 
bility distribution over |p PT |, where p is a generic real two-qubit density matrix. The domain of 
separability is \p PT \ > 

V. CONCLUDING REMARKS 



It is clear that it would be of considerable utility to have available exact values for 
still higher-order (than m — 9) moments-and for the coefficients C2j(m) of the terms in 
the intermediate functions J m (/x)-but the associated computational demands seem quite 
considerable. (The Hilbert-Schmidt separability probability predicted by the curve in Fig. 0]- 
le "probability mass" lying within the interval [0, o|g]— is 0.39648, while our previous 



that is t 
studies 



19j , indicate that the actual value is considerably higher, « 0.45-a discrepancy still 



higher-order moments should ameliorate.) Also, of course, it would be interesting to extend 
our forms of analyses from the real case to the more fully general setting of the 15-dimensional 
convex set of complex two-qubit (4 x 4) density matrices. But, at this stage of development 
of our technical apparatus, we are unable even to compute the corresponding first Hilbert- 
Schmidt moment (mean) over \p PT \ (known to be — ^ in the real two-qubit case). However, 
progress, in these regards, should be much more readily achievable in terms of the moments 
over \p\. 

In fact, inputting the first twenty such moments computed into the FindSequenceFunction 
command of Mathematica, we obtained the formula (again employing the Pochhammer 
symbol) 

256~ m (l) m (2) m (3) m 108972864000r(m + l)Y{m + 2)T(m + 3)T(m + 4) 



c 



m/ 'complex 



(ii) (2) (m 

\4/m\2/mUit 



T(4m + 16) 



(58) 
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The associated moment- generating function is the generalized hypergeometric function 

M(t) complex = 3 F 3 ( 1,2,3; — , -, — ; — J . (59) 

Similarly, for the quaternionic two-qubit case, 

~, 315071454005160652800000r(m+ l)r(m + 3)r(m + 5)r(m + 7) 

(m/quat ~ I\4(m + 7)) ' ^ ' 

and 

/ 29 1 5 / \ 

M(t) quat = 3 F 3 1, 3, 5; -f , ^, ^; -M (61) 



4 ' 2 ' 4 ' 256, 

We have also been able to compute exactly the first twelve moments of the probability 
distribution of |p| over [0,6~ 6 ], where p is a generic complex qubit-qutrit (6 x 6) density 
matrix, but this seemed to be an insufficient number of moments to discern a general formula. 
All eleven moments found were the reciprocals of positive integers. The first moment (mean) 
was 449 g 388 = (2 2 • 3 ■ 13 • 19 • 37 ■ 41) _1 , while the second moment was smaller than this by a 
factor of 153 g 939 . (The twelfth moment is approximately 6.16876 • 10 ~ 68 .) 
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